Analytical solutions for the dynamics of two trapped interacting ultracold atoms 
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We discuss exact solutions of the Schrodinger equation for the system of two ultracold atoms con- 
fined in an axially symmetric harmonic potential. We investigate different geometries of the trapping 
potential, in particular we study the properties of eigenenergies and eigenfunctions for quasi-one- and 
quasi-two-dimensional traps. We show that the quasi-one- and the quasi-two-dimensional regimes 
for two atoms can be already realized in the traps with moderately large (or small) ratios of the 
trapping frequencies in the axial and the transverse directions. Finally, we apply our theory to 
Feshbach resonances for trapped atoms. Introducing in our description an energy-dependent scat- 
tering length we calculate analytically the eigenenergies for two trapped atoms in the presence of a 
Feshbach resonance. 
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I. INTRODUCTION 

Atomic interactions represent one of the major ingredi- 
ents for the schemes implementing quantum information 
processing in systems of neutral trapped atoms. Devel- 
opment of optical lattice technology [l|, atom chips 
and dipole traps [3j allows to create tight external con- 
finement for neutral atoms. Feshbach resonances, widely 
used in recent experiments on ultracold atoms, permit 
for tuning of atomic interactions, which was the key in- 
gredient to accomplish molecular Bose-Einstein conden- 
sates and the superfluidity of fermionic atoms in ultra- 
cold gases 4] . In addition, realization of Mott insulators 
[3 gives the possibility to precisely control a number of 
atoms confined in a single well. All these achievements 
make systems of ultracold neutral atoms very attrac- 
tive in the context of quantum information processing or 
quantum control at the atomic level. Moreover, they also 
open a way for experimental studies of few-body interact- 
ing systems in the presence of tight external potentials. 

On the theoretical level, the system of two interact- 
ing atoms in a harmonic trap can be solved analytically 
for spherically symmetric [|| or axially symmetric har- 
monic potentials Q- Both these approaches model the 
interaction in terms of an s-wave delta pseudopotential 
H, 13 • Generalization of the pseudopotential to higher 
partial waves \WL UA\ allows to solve the problem for 
generic types of (central') interactions, in the presence of 
spherically symmetric J10l | or axially symmetric [ljj har- 
monic traps. Moreover, supplementing the pseudopoten- 
tial with an energy-dependent scattering length jl^, uM 
extends the validity of the analytic results to the case of 
very tight traps or large scattering lengths, and accounts 
properly for the whole molecular spectrum [l(j- Such 
model provides for a very accurate description, which has 
been verified, for instance, in the recent experiment on 
the creation of molecules of fermionic atoms in an optical 
lattice [bH- 



In this paper we discuss in detail the exact solutions 
for two interacting atoms confined in axially symmetric 
harmonic traps Jg . We present derivations of the ana- 
lytical results discussing different geometries of the trap- 
ping potential. In the limiting cases of the quasi-one- 
and quasi-two-dimensional traps the system can be ef- 
fectively described in terms of a lower-dimensional trap 
with renormalized scattering length. We investigate the 
limits of applicability of the quasi-one- and the quasi- 
two-dimensional descriptions showing that they are valid 
already for moderately large (or small) ratios of the trap- 
ping frequencies in the axial and transverse directions. 
Finally, we consider the effects of Feshbach resonances 
on the trapped atoms. Employing the standard theory 
of Feshbach resonances we express the energy-dependent 
scattering length in terms of the usual parameters de- 
scribing the resonance. This allows to derive an explicit 
formula determining the energy spectrum in the presence 
of a Feshbach resonance. 

The paper is organized as follows. In section ^ we 
derive analytical solutions of the Schrodinger equation 
for two atoms confined in an axially symmetric trap, 
by expanding them in the basis of harmonic oscillator 
wave functions. Section Hill is devoted to the analysis 
of eigenenergies. In particular, in section Hit Al we derive 
analytical results for cigar-shape traps, while the limiting 
case of the quasi-one-dimensional traps is analyzed in sec- 
tion llll Bl Section llll Cl presents the analytical results for 
pancake-shape traps. The quasi- two-dimensional regime 
in these traps is studied in section llll Dl Section llVl ana- 
lyzes the properties of eigenfunctions. In section llV Al we 
derive two series representations for the wave functions, 
which are valid for arbitrary ratio of radial to axial trap- 
ping frequencies. The behavior of the eigenfunctions in 
quasi-one- and quasi-two-dimensional traps is discussed 
in sections llV Bl and IIV CI respectively. In section Ivl we 
illustrate the applicability of our theory, calculating the 
energy spectrum for two atoms interacting in the pres- 



ence of a Feshbach resonance. We end in section IVTI 
presenting some conclusions. Appendix [X] presents some 
technical details related with the derivation of the energy 
spectrum for quasi-two-dimensional traps. 



II. SYSTEM 

We consider two interacting atoms of mass m confined 
in an axially symmetric harmonic trap with frequencies 
10 z and u>± in the axial and transverse directions, respec- 
tively. The Hamiltonian of the system reads 



H 



2m 



v?- 



2m 



V| + K(ri) + l4(r 2 ) + M(r 1 -r 2 ), (1) 



where ri and r 2 denote the positions of the two atoms, 
and T4(r) is the trapping potential 



V t (r) = —(u ± p 2 +u z z 2 ), 



(2) 



with p 2 



V 



For sufficiently low energies, the 



scattering is purely of s-wave type and we model_the 
interaction potential by a regularized delta function 



4:Trh 2 a d 

Vi{v) = S(r)—r, 

m or 



(3) 



with a denoting the s-wave scattering length. 

For the harmonic trapping potential, the center-of- 
mass and relative motions are decoupled. Substituting 
new coordinates r = ri — r 2 and R = (r x + r 2 )/2, we 
decompose the total Hamiltonian into the center-of-mass 
part Hqm and relative part -ff re i: 



Hqm 



-ffrel = 



2M 

h 2 



V 2 



M 



Vt(R) 



Ip, m 



(4) 
(5) 



where p = m/2 and M = 2m denote the reduced and 
the total mass, respectively. 

In the following we use dimensionless variables, in 
which all lengths are expressed in units of d = y/h/(/j,uj z ), 
and all energies are expressed in units of fttu z . The eigen- 
functions of the center-of-mass Hamiltonian H cm a re the 
usual harmonic-oscillator wave functions. The eigenfunc- 
tions \P(r) of the relative motion have to be determined 
from 



1 r 
-2^ 



1 d ' 

-(r] 2 p 2 + z 2 )+2iraS{r)—r 



*(r) = £*(r) 

(6) 

where r\ = ujj_/uj z . Its solutions can be found by de- 
composing \&(r) into the complete set of the harmonic 
oscillator wave functions |6( 



*(r) =^2cn t k$ n ,o(p,v)9k(z) 



(7) 



where, $ n ,m(p, f) denotes the states of the two- 
dimensional harmonic oscillator in polar coordinates 
(p, (p) with the radial and angular quantum numbers 
n and m, respectively, whereas 0&(z) is the one- 
dimensional harmonic-oscillator wave function with the 
quantum number k. We note, that only the states with 
m = 0, i.e. the states with vanishing angular momentum 
along z enter the summation in Eq. (JJJ . The states with 
m 7^ vanish at r = 0, and they are not perturbed in the 
presence of the interaction potential. Substituting the 
expansion J7J) into the Schrodinger equation © yields 

= Y^ Cn, k (E nik - E)® n , {p, v)Qk(z) (8) 

a . k 

d 

+ 2TraS{r)—r^2c n ^ ny0 (p,if)e k (z), (9) 



where E n ^ = l/2+rj + k+2r]n are dimensionless eigenen- 
ergies of the three-dimensional axially-symmetric har- 
monic oscillator. To determine the expansion coefficients 
c n,k we project Eq. © onto state $ n ' o(p, <p)&k'(z) with 
arbitrary nl and k', obtaining 



^n,fe ^ 



*n,o(o.^)©3E(o) 



En 



E 



(10) 



where C is a constant fixed by the normalization of the 
wave function. The value of C is related to the expansion 
coefficients c„ & through 



C = 2na 



8 



o- ry^Cn,fc$ w ,o(p, <p)9k{z) 



Or 



(11) 



Substituting the solution I jlOj l for coefficients c n ^k into 
Eq. 1|11[1. the numerical constant C disappears, and we 
obtain an equation which determines the eigenenergies 
with m = 0: 



1 
2na 



d_ 

Or 



r* £ (r) 



(12) 



r=0 



where 



*fW = 2^— n„„ I I c (I 3 ) 



n ,k 



2rjn + k — £ 



n , k 



and £ = E — Eq denotes the energy shifted by the zero- 
point oscillation energy Eq = 1/2 + r\. For values of £ 
solving Eq. i|12|l . the functions \&£(r) represent the non- 
normalized eigenstates of the relative motion. We note 
that \&£ (r) is proportional to the single-particle Green 
function G(0, r) of the three-dimensional anisotropic har- 
monic oscillator. We stress that, the regularization oper- 
ator -j^r and the summation in Eqs. I|12l) and (|13fl cannot 
be interchanged, and the summation over n and k must 
be done first. This is related with the divergent behavior 
of the Green function at small r, which is regularized by 



■gpr. To perform the summation in Eq. i|13|) we express 
the denominator in Eq. (|13fl m terms of the following 
integral: 



2rjn + k — £ 



d< e ~t(2 V n+k-£) ^ 



(14) 



Since n > and k > 0, the integral representation 
(|14H is valid for £ < 0. The wave function of the two- 
dimensional harmonic-oscillator with m — is given by 



*„,0(P,V) = ^e-i>> a /*L n (T,f?) 



(15) 



where L n (x) is the Laguerre polynomial. To perform 
the summation over n in Eq. (j!3|> . we utilize the fact 
that L n (0) — 1 and we apply the generating function for 
Laguerre polynomials |l6| 



Em 



X)Z 



(1 — z) 1 exp 



1 



(16) 



On the other hand, the sum over k involves the eigen- 
states of the one-dimensional harmonic-oscillator 



e fc (2) 



-z 2 /2 



7rV*v^Jbi 



#*(*), 



(17) 



where Hk{z) is the Hermite polynomial. The summation 
can be performed with the help of the following gener- 
ating function for the products of Hermite polynomials 



na 



fc=0 



2 k kV 



E^r^)^) 



Vl~¥ 



(18) 



Upon inserting Eq. I|14|) into (|13|l . and performing the 
summation according to (|16fl and (|18|) . we find the fol- 
lowing integral representation for *&£ (r) 



*£(r) 



(27T)* 



df 



exp 



t£ 



• coth t 



2§- coth(7?t) 



•^/sinh(i) sinh(?7t) 



(19) 



The integral is convergent for E < Eq, however, the valid- 
ity of the final result will be extended to energies E > Eq, 
by virtue of the analytic continuation. 

Let us investigate now the behavior of '$>£ (r) for small 
values of r. In the limit r — ► 0, the main contribution to 
the integral comes from small arguments t. In the leading 
order we can neglect the dependence on the energy E, 
and the expansion of Ijl9|l for small t yields 



*f(r) 



1 



(27T) 



3/2 



dt 



-r 2 /(2t) 
t 3 / 2 



1 

2nr ' 



r< 1. 



(20) 

We note, that, for small r the function ^(r) diverges in 
the same way as the Green function in a homogeneous 



space. This is related with the fact that at short distances 
the behavior of the wave function is determined mainly 
by the interaction between the particles. After extracting 
the divergent behavior of ^>£ (r) , we can simplify Eq. i|12|) 
determining the cigcncncrgics. To this end, we substitute 
the integral representation l|19|l into (|12|) . and subtract 
from $>£ (r) the r.h.s. of H2U|) . This can be done since the 
term l/(27ir) is removed by the regularization operator, 
and does not give any contribution to JT2J. In this way 
we obtain a simpler expression determining the energy 
levels for m = 



a 



where 



T(x) 



At 



r\e 



VI - e~t (1 - e- 



nt\ 



1 

73/2 



(21) 



for x > 



(22) 

Similarly to l|19|) , the validity of the integral representa- 
tion IL'l't is limited to £ < 0. In general, T(x) can be 
calculated from the following series representation, valid 
for all values of x: 



F(x) 



n 



E 



T(x + nrj) 



1 



27T ^oV^ + x + n v) VvV^+T 



Vv 

2vr 



c(|), 

(23) 
where T(x) is the Gamma function and C(i) denotes the 
Riemann zeta function. To derive Il2.''>[| we retain the 
summation over n in Eq. I|13|l. while we perform the in- 
tegration over t, and the rest of the steps is the same as 
in derivation of (|22|1 . 

Finally, in the numerical calculations of !F(x) we were 
using a recurrence relation, which relates the values of 
function T(x) for different arguments x. To derive this 
formula we calculate the difference between J-(x) and 
T[x + rj) obtaining with the help of i|22|) 



f°° ne~ xt T( 

T{x) - F(x + rj) = / dt ' = vV^t^ 

Jo vl-e - * r(x- 



T(x)_ 

2> 



(24) 



III. EIGENENERGIES 

A. Cigar shape traps 

In this section we discuss the energy spectrum of two 
interacting atoms confined in a trap with rj > 1. Let us 
first consider the case r\ = n, with n being a positive in- 
teger. In this particular case we can treat the function 
1 — e~ nt , which appears in the denominator of the in- 
tegrand in Eq. (|22|l . as a polynomial of order n in the 
variable e _t , and perform a decomposition into a sum of 
simple fractions 



_. n— 1 



1-e- 



1 ( D — t—i27rmfn 



(25) 



Substituting this decomposition into 1(22(1 yields 



T(x) 



OO 

df 



n _ e -*) 3/2 t 3/2 

/"OO 



n-1 



,-i •'0 



-ait 



VI - e~* (1 - e-'- 



■i27rm/nA 



= -2V^ 



r(x) 



r(x-|) 



ra-l 



v 2 7 m— 1 



.1. ..<.■.. J | r) 5 C 



(26) 



where F (a, 6; c; x) denotes the hypergeometric function. 
We note that Eq. 1(26(1 is derived from the integral repre- 
sentation 11' I'D applicable for x > 0, however the validity 
of the final result can be extended to x < by virtue 
of the analytic continuation. Despite the presence of the 
complex roots of 1 in the argument of the hypergeometric 
function, it can be easily verified that the whole expres- 
sion remains real for x € R. For the special case n = 1, 
the second term in Eq. (|26|l disappears and we obtain 
the well-known result for the spherically symmetric trap 

una 



F{x) = -lypH 



r(s) 



n 



(27) 



Fig. ^ presents the energy levels calculated for 77 = 5 
from Eq. 121(1 . with F(x) given by the exact formula 
(126(1 . For a — the eigenvalues are given by the poles of 
.F(x), and we recover obviously the energy spectrum of 
the harmonic oscillator. On the other hand, the eigen- 
values for a — ► +00 and for a — > — 00 approach the same 
asymptotic values, corresponding to zeros of F(x). The 
level spacing for large a is not uniform, as in the case of a 
spherically symmetric trap |6j , and the distance between 
energy levels is larger every fifth level, which results from 
the geometry of the trap. For a > (repulsive poten- 
tial) the energy levels are shifted upward with respect to 
the non-interacting case, while for a < (attractive po- 
tential) they are shifted downwards. In a homogeneous 
space, the three-dimensional regularized delta potential 
possesses a single bound state for a > 0, with energy 
E = —h 2 /ma 2 . From Fig. ^ we see that such a state is 
also present in the case of harmonic confinement, how- 
ever, its energy is shifted upward due to the presence 
of the trap. Moreover, the presence of external confine- 
ment gives rise to the appearance of a bound-state also 
for a < 0, which can be observed in Fig. Q]as a branch 
of the spectrum starting from the energy of zero-point 
oscillations. 

Fig. [2] presents the exact energy levels for 77 = 1.0 (up- 
per plot) and 77 = 1.1 (lower plot). The former are given 
by the analytical result 1(27(1 . while the latter were cal- 
culated numerically from Eq. (1221) for £ < 0, and with 
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FIG. 1: Eigenenergies of the relative motion for two atoms 
interacting via s-wave pseudopotential and confined in a har- 
monic trap with 77 — lu±/uj z = 5. The scattering length 
a is scaled in units of the harmonic oscillator length d = 

s/h/({iuj z ). 



the help of the recurrence relation ((24(1 for £ > 0. Com- 
paring upper and lower plots, we note that the energy 
spectrum for 77 = 1.1 has a richer structure than the one 
for the spherically symmetric trap (77 = 1). In the latter 
case some of the excited states are degenerate, and they 
do not appear in the energy spectrum given by Eq. 1J2*T|| . 
As can be easily verified by taking linear combinations of 
degenerate wave functions, the number of solutions not 
vanishing at r = can be always reduced to one, and as a 
consequence only one of the degenerate states is affected 
by zero-range interaction. Analyzing the behavior of the 
eigenenergies close to a = 0, we notice the appearance of 
avoided crossings for non-integer values of 77. 



B. Quasi- ID regime 

Now, we turn to the discussion of the energy spectrum 
for 77 3> 1 . For x ~ 77 and 77 3> 1 , the main contribution to 
the integral in 1(2211 comes from small arguments t. In this 
regime we perform the approximation y/l — e~ l « \ft in 
the denominator of ((22(1 and obtain 



T(x) 



OO 



■qe 



Vt(i 



-vt) 



1 

W 1 



, 7? »1 (28) 



Now the integration can be done analytically. This yields 



F(x) w 0t77Ch(1/2,x/?7), x>?7 



(29) 



where Cn{s,a) denotes the Hurwitz Zeta function: 
( H (s,a) = J2kLo( k + a )~ S 2l V To extend the valid- 
ity of the approximate result ((29(1 to x positive and much 
smaller than 77, or to negative x, we make use of the re- 
currence relation l(24|) . Applying the recurrence formula 
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FIG. 2: Eigenenergies of the relative motion for two atoms 
interacting via s-wave pseudopotential and confined in a har- 
monic trap with r\ — uj±/uj z = 1.1 (upper plot) and r) — 1 
(lower plot) . The scattering length a is scaled in units of the 
harmonic oscillator length d — y/h/{^iw z ). 



once, we find 



T{x) w 0F7/ C,h ( |,1 



A r- r («) 

W ' r(x + 



(30) 



Numerical comparison of the exact result l|2tj[l with the 
approximation (|3U[) shows that the latter provides quite 
accurate values for x > —r\. Thus, utilizing i|30|l we are 
able to calculate the energy of the ground state and of 
the first excited states, up to £ < 2r\. To find the energies 
of higher excited states one can apply recursively H24|) . 

It is interesting to compare our results with the pre- 
dictions of one-dimensional model, where the scattering 
length is renormalized due to the tight confinement in the 
transverse direction 18]. At low energies of the scattered 
particles, the one-dimensional pseudopotential takes the 
form 18] 



V 1D (r) 



MOID 



-S(r) 



(31) 



where am is the one-dimensional scattering length. Re- 
peating similar steps as for the three-dimensional trap, 
we obtain the following implicit equation determining the 
eigenenergies of the relative motion in a one-dimensional 
harmonic potential 6] 



2aiD = 



r(- 



r(- 



(32) 



We compare the latter equation with the three- 
dimensional result for r/ 3> 1, 



1 



= 07 Ch 



1- 



£_ 

2 V 



i] 



r(-f) 
■! + *)' 



r(- 



(33) 



obtained by substituting into Eq. (|21(l the function T(x) 
given by l|3U|l . For the energies \£\ <C r\ one can ne- 
glect the dependence on energy in the first term on the 
r.h.s. of (El: Ch (1/2, 1 - £/(2??)) « C(V 2 )- Then it is 
straightforward to observe that Eqs. l|32"|l and l|35|) give 
the same energy spectrum, provided the one-dimensional 
scattering length am is related to a by 



am = 



1 
2r)a 



C(i) 



(34) 



Expressing this relation in physical units, we recover the 
result of Ref. [3 



am 



2a 2 ^ 2J ' 



(35) 



where d± = ^/h/(puj±). 

In the previous works [3, Il5j it was shown that a one- 
dimensional model with renormalized scattering length 
provides a very accurate description of the spectrum for 
E > E Q . On the other hand, this approach is not valid 
for energies E < E when the system possesses a bound 
state. Its energy has to be determined from Eqs. I|21|) 
and (|29|l . derived from the three-dimensional approach. 
The condition for the energy of a bound state expressed 
in the physical units reads 



A = c (- E °~ E 

a \2' 2huj± 



(36) 



which agrees with the results of [ijj- We note that 
Eq. H36|) involves only the trapping frequency iu± in the 
tightly confined direction. Due to this reason the one- 
dimensional model, which depends crucially on w Z} fails 
to describe the energy of a bound state. 

We have verified that the quasi-one-dimensional regime 
for two interacting atoms does not require r\ very large, 
but is already realized for ij ~ 10. The case of -q = 10 
is illustrated in Fig. |3| where we compare the exact en- 
ergy levels with the ones calculated from Eq. (|32*|l for 
the one-dimensional spectrum with <2id given by H34fl . 
We observe quite good agreement for the lowest eigen- 
states. For higher excited states the two approaches 
start to differ around the unitarity point (1/a = 0). 
The predictions of the one-dimensional model can be im- 
proved when the eigenenergies are calculated assuming 
an energy-dependent one-dimensional scattering length 

m 



fm d\ d ± (I Eo 



E 



(37) 



This expression can be obtained from the result l|29|l : 
however, in this case we do not apply the approximation 
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FIG. 3: Eigenenergies of the relative motion for two atoms 
interacting via s-wave pseudopotential and confined in a har- 
monic trap with r\ = lj±/uj z — 10. The exact energy lev- 
els (solid lines) are compared with predictions for the one- 
dimensional model with energy-dependent (dotted line, in- 
distinguishable from the solid one) and with the standard, 
energy-independent renormalized scattering length (dashed 
lines). The three-dimensional scattering length a is scaled 
by the harmonic oscillator length d = ^Jh/(iMJ z ). 



£ = for the first term of l|29|) . We note an excellent 
agreement between the exact energy levels and the one- 
dimensional spectrum evaluated with an)(E), which on 
the scale of Fig. are indistinguishable. 



C. Pancake shape traps 

In this section we investigate the energy spectrum of 
the two interacting atoms confined in harmonic traps 
with < r] < 1. First, we derive an explicit formula 
for ^(x) in the case when r\ is the inverse of a positive 
integer. We start from the integral representation H22|) 
of function T(x), which for r\ = 1/n takes the following 
form: 



T{x) = / dt 
Jo 



i vT 



'i 



,-t/i 



) 



l 
IN* 



(38) 



To calculate the integral we make use of the following 
identity 



1 



,—t/r, 



1 n — 1 

1 y 

l-e-* ^ 



-tm/r, 



(39) 



To prove this identity one can use the formula 1 



(1 - ie)(1 + x + x 2 + ... + X 11 - 1 ) with x 
tuting ||2H1) into JSSJ leads to 



"*/". Substi- 



F{x) 



1 71 — 1 



dt 



-t(x-\-7n/n) 



1 
^372 



(40) 
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FIG. 4: Eigenenergies of the relative motion for two atoms 
interacting via s-wave pseudopotential and confined in a har- 
monic trap with 77 = uj±/uj z — 1/5. The scattering length 
a is scaled in units of the harmonic oscillator length d — 
s/K/(jmJz). 



which can be evaluated analytically 



T(x) 



20F" 



1 



n f 



to n* - 



7X1 \ 

a ' 



(41) 



We note that for the special case of a spherically sym- 
metric trap (n = 1), we recover obviously the result l|27|l . 
The exact energy spectrum calculated from combined 
Eqs. (|2"TJl and l|H")l. and for 77 = 1/5 is presented in Fig.0] 
We observe that the energy spectrum has similar features 
as for the cigar shape traps. In the limit of a — > ±00 the 
eigenenergies approach the same asymptotic values, irre- 
spective of the sign of the scattering length. In contrast 
to the cigar shape traps, the level spacing at a — > ±00 is 
almost uniform, except for the gap between the ground 
and the first-excited state. Again, due to the external 
confinement, the solutions with energy E < Eq, corre- 
sponding to the bound states of the interaction poten- 
tial, occur both for positive and negative values of the 
scattering length. 



D. Quasi-2D regime 

Now, we turn to the analysis of the energy spectrum in 
the quasi- two-dimensional traps with r\ <C 1. For simplic- 
ity we assume r\ = 1/n with n being an integer, however, 
this choice is not crucial for the applicability of the final 
results, as we show later. We start our derivation from 
formula (|41|l . in which for n » 1 we approximate the 
summation by an integration 
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[dtB(a 
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(42) 



where B(x, y) = T(x)T(y)/T(x + y) is the Euler beta 
function. This approximation is valid for x > 0, which 



guarantees that the function B (x + t, —5) is free from 
the singularities in the interval of integration. As we 
show in Appendix EI the integral in Eq. (|43|) can be 
expressed in terms of the following functions: 



I At B(x + t,-\) = -$(z) -\\\x, 
Jo 



(43) 



with 
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(44) 



Eq. 12111 together with the approximate result l|42[) de- 
termine the energy of the bound state {£ < 0) in the 
quasi-two-dimensional traps. Expressing this result in 
the physical units we obtain 
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We observe that in quasi-two-dimensional traps the prop- 
erties of a bound-state depend solely on the trap fre- 
quency in the tightly confined axial direction. For a shal- 
low bound state (Eo — E <C hu> z ) we can approximate 
$(— £/2) by $(0) w 1.938, and in this regime we recover 
the result of Ref. 0: E Q - E = hu z 0.288exp{^d/a). 
Let us investigate now the energy spectrum for £ > 
containing the excited states. In the following, we will 
compare our results, obtained from the three-dimensional 
description with predictions for the two-dimensional sys- 
tem, with the scattering length renormalizcd due to the 
tight confinement in the axial direction [221 123| . In anal- 
ogy to the derivation of Huang and Yang [ij, one can 
show that in two dimensions the s-wave pseudopotential 
takes the form 

V 2D (r) = r^J r*(r) ( 1 - HAkp)p^-) , (46) 

/iln(fca 2 D) V UP J 

where k 2 = 2pE/h 2 , A = e 7 /2, with 7 denoting the 
Euler constant and a 2 D is a two-dimensional scattering 
length related with the the s-wave scattering phase shift 
So by tan<5o = (tt/2) ln~ {ka 2 o)- The regularization op- 
erator (1 — In(Akp)pJp) removes the logarithmic-type di- 
vergence of the two-dimensional scattering solution at 
p = 24]. We note that even in the limit k — > 0, the 
pseudopotential (|4^|l depends on energy. 

Using similar techniques as for the three-dimensional 
trap, one can show that the energy spectrum of two inter- 
acting atoms confined in the two-dimensional harmonic 
trap of frequency uj± is given by |(j 
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(47) 



where ip(z) denotes the digamma function: ip(z) — 
(d/dz) lnT(z). Expressing Eq. I|47|l in terms of dimen- 
sionless units related with w z , we obtain 
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(48) 



To find the connection between the two-dimensional and 
the quasi-two-dimensional energy spectrum we use the 
following approximate formula, valid for fjCl (see Ap- 
pendix El for derivation) : 

F(x) w -$(«) -hiTj- i){x/rj), 77 -Cl (49) 

Substituting (|49|l into J2U, we obtain the condition 
for eigenenergies of excited states in the quasi-two- 
dimensional traps: 

&=]n(r,) + 9(-£/2)+1> (-£/&)) (50) 
a 

For the lowest excited states (\£\ <C 1), we can neglect 
the energy dependence of $(— £/2). In this regime it 
is straightforward to observe that Eqs. (|15|l and l|5*U|) 
predict the identical energy spectrum, provided that the 
two-dimensional scattering length a 2 u is related with a 
by 
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" 2 " = 7I CXP v— ~^ 



(51) 



Consequently, the two-dimensional coupling constant 
<72D = — ttH 2 / (p\ji(ka2i))) expressed in the physical units 
is given by 
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(52) 



which agrees |2jj with the results of Ref. |22ll23|. 

Similarly as for the quasi-one-dimcnsional traps, we 
have also verified the limits of applicability of the two- 
dimensional model with renormalized scattering length. 
It turns out that the latter approach is quite accurate up 
to 77 ~ 1/10. We note that in the case of r\ = 1/10 the ra- 
tio between the harmonic oscillator lengths in the radial 
and axial direction is equal to vlO, therefore we would 
expect that the shape of the wave function is rather far 
from the quasi-two-dimensional one. Hence, the good 
agreement between the exact energy levels and the two- 
dimensional spectrum is quite surprising. This feature is 
shown in Fig.|3 We observe that the two-dimensional ap- 
proximation gives less accurate results for higher excited 
states around the unitarity point (1/a = 0). This again 
can be improved by introducing an energy-dependent 
two-dimensional scattering length 
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which can be easily derived by comparing i|50[) with 14811 . 
Eigenenergies calculated assuming the energy-dependent 
scattering length o, 2 d(E) are in excellent agreement with 
exact ones, which is illustrated in Fig. [S] 




FIG. 5: Eigenenergies of the relative motion for two atoms 
interacting via s-wave pseudopotential and confined in a har- 
monic trap with r\ — uj±/uj z = 0.1. The exact energy lev- 
els (solid lines) are compared with predictions of the one- 
dimensional model with the energy-dependent (dotted line, 
indistinguishable from the solid one) and with the standard, 
energy-independent renormalized scattering length (dashed 
lines). The three-dimensional scattering length a is scaled 
by the harmonic oscillator length d — tJH/((au) z )- 



IV. WAVE FUNCTIONS 



A. Axially symmetric trap of arbitrary anisotropy 



Let us turn now to the analysis of the wave functions. 
For m y^ the wave functions of noninteracting particles 
vanish at r = 0, and as a consequence they are not mod- 
ified by the presence of the zero-range potential. The 
nontrivial m = wave functions are given by Eq. (|13H . 
or by the integral representation (|19fl . valid for £ < 0, 
where in the place of £ one has to substitute eigenener- 
gies calculated from Eq. (|21|l . In the following we derive 
analytic formulas for the wave functions which are sim- 
pler and more convenient for numerical calculations. 



We start from Eq. (|T5j) : 
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. (54) 
where we have substituted explicitly harmonic-oscillator 
wave functions. In the next step we rewrite the factor 
(2i]m + k — £) in the denominator of (|54|l as the integral 
HI 41 . Then, we perform the summation only over a single 
variable k or m, using the generating functions l|16fl and 
(|18fl . This yields an integral, which can be calculated 
analytically, and the final result can be expressed in terms 
of a series in a single variable. 

Let us first perform the summation over the quantum 



number k. Applying the generating function (|18|l we get 
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(55) 



The integration can be done analytically by introducing 
a new variable of integration x = e~ 2t /(l — er 2t ) 16]. In 
this way we find 
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(56) 



m=0 



where D v {x) denotes the parabolic cylinder function. 
Since the parabolic cylinder functions D v (x) are well 
denned both for positive and negative values of index 
u, the validity of the result (|56(l is automatically ex- 
tended to all values of energy £. We note that parabolic 
cylinder functions -D„(|z|v2) are one-dimensional wave 
functions for two interacting atoms in a harmonic trap. 
In the transverse direction the expansion l|56l) involves 
harmonic-oscillator wave functions. The latter constitute 
an orthonormal basis, which simplifies analytic calcula- 
tions of the matrix elements involving the wave functions 
in position representation. As an example let us calculate 
the normalization factor Af for the wave function ^(r) : 
N~ 2 = Jd 3 r|* £ (r)| 2 , where for ^s(r) we substitute the 
expansion (|56|l . Integration over p is trivial due to the 
orthonormal properties of the transverse wave functions, 
whereas the axial integration involving D v (\z\\/2) 2 can 
be performed analytically. This results in 
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with (3(x) = [ip((x + l)/2) - ip(x/2)]/2. 

Another expansion is obtained when in Eq. (|54|) we 
perform a summation over the quantum number m. Rep- 
resenting the denominator of H54|) as the integral (|14|) and 
utilizing the summation formula (|16f) . we arrive at 
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The result of the integration can be expressed in terms 
of the confluent hypergeometric function f/(a, 6, z), which 



yields 
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(59) 



where we have substituted the values of the Hermite 
polynomials at z = 0: H 2k (0) = (-2) fc (2fc - 1)!!, and 
-ff2fe+i(0) = 0. As in the previous case, this expansion 
can be utilized both for positive and negative values of 
energy £, since the analytic continuation is provided au- 
tomatically by the properties of T(a) and U(a,l,x) for 
negative values of the parameter a. As it can be eas- 
ily verified, U(a, l,rjp 2 ) is a wave function of the two 
interacting atoms in a two dimensional trap. In the ax- 
ial directions, the expansion is done in the basis of one- 
dimensional harmonic oscillator wave functions. We note 
that U(a, 1, rjp 2 ) exhibits logarithmic divergence at p = 
and the result l|59[l is not correct for p = 0. In this par- 
ticular case the expansion ij5*6"jl should be used instead. 

The expansion (|59|) can be used to obtain the nor- 
malization factor JV. Substituting l|59|) into JV~ 2 = 
J d 3 r|\E r £(r)| 2 and performing first the integration in the 
axial direction, and then applying 



e- x dxT{a) 2 U(a,l,x) 2 =£ H (2,a) (60) 



for the transverse integration we arrive at 
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In several applications of our theory, like evaluation of 
the exact dynamics for two interacting atoms, it is nec- 
essary to determine the complete set of eigenfunctions. 
The eigenfunctions that are modified by the interaction 
are given by representations (|56|l and l|59|) . with the en- 
ergy £ evaluated from the implicit equation 121|) . The re- 
maining eigenfunctions are the harmonic-oscillator wave 
functions that vanish at r = 0. Applying the notation of 
Section ITU thev can be written as $> n ,m(p, <p)®k{z) with 
m > or m = and k odd. When uo z and uj± are in- 
commensurate all the rest of the eigenstates is generated 
by (|21|l . In the case of commensurable trapping frequen- 
cies, like for example in the particular case of 77 = n or 
j] = 1/n, there can appear accidental degeneracies in the 
energy spectrum of the harmonic oscillator. Assuming 
that there are N degenerate eigenstates, the condition 
(|21|l generates only the single state affected by the in- 
teraction, while the remaining N — 1 states can be de- 
termined, for instance from the Gram-Schmidt orthonor- 
malization procedure, and they will automatically vanish 
at r = 0. 



B. Quasi-ID regime 

Now we discuss the behavior of the wave functions in 
the limit of very elongated traps: 77 3> 1. We start from 



the integral representation i|19|) to obtain the wave func- 
tion of the ground state. For 77 ^> 1 and \£\ ~ r\ the 
main contribution to the integral comes from the region 
of small t. Expansion of the integrand to the lowest order 
in t yields 
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Changing the variable of integration and expressing the 
final result in the physical units we obtain 
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We note that the wave function of a bound state depends 
only on the trapping frequency lo± in the tightly confined 
direction [2!j . This feature has been already observed on 
the level of eigenenergies. 

The wave function of the bound state in quasi-one- 
dimensional traps can be also calculated from the follow- 
ing expression: 
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(64) 

To derive (|64|l we expand l|55l) for small t and perform 
analytically the integration over t. From Eq. I|64l) one 
can easily derive approximate axial and radial profiles of 
the bound-state wave function 0] , which agrees very well 
with the exact one evaluated from Eqs. I|56|l and l|59|l . 

The exact ground-state wave function for a = ±00 
and 77 = 100 is presented in Fig. H3 For small r it is al- 
most isotropic, due to the presence of the divergent factor 
l/(27rr), while for larger values of r the wave function is 
slightly elongated in the z-direction, which reflects the 
geometry of the trap. 

Let us discuss now the properties of the excited-state 
wave functions. In the regime of energies corresponding 
to the lowest excited states and for the arguments z not 
too small, it turns out that the first term of the series H56|) 
dominates the sum. This is a consequence of the behavior 
of the parabolic cylinder function D v (x), which are fast 
decaying when index v becomes large and negative. Thus 
the approximate wave function of the excited state reads 
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Obviously, the latter approximation is not valid for small 
r, where the exact wave function exhibits a divergent be- 
havior: 1/(27T7"). Fortunately, the main part of the wave 
function is located for larger values of r, and the region 
of small r gives a rather small contribution to the total 
wave function. On the other hand, we should keep in 
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FIG. 7: Exact wave function r*$!(r) for two atoms interacting 
via s-wave pseudopotential and trapped in a harmonic po- 
tential with rj — w±/u) z = 100. The figure presents the first 
excited state for the scattering length a — ±co. All lengths 
are scaled in units of a z = \/h/ 'puj z . 



FIG. 6: Exact wave function r*l/(r) for two atoms interacting 
via s-wave pseudopotential and trapped in a harmonic po- 
tential with rj = uj±/uj z — 100. The figure shows the ground 
state for the scatter ing lengt h a — ±oo. All lengths are scaled 
the units of a z = ^/h/{pu z ). 



mind that the delta pseudopotential is only an approxi- 
mation, and for small r, comparable to the effective range 
of the physical potential, the behavior of the real wave 
function is quite different from the predictions based on 
the pseudopotential. 

To estimate the accuracy of the approximation i|65|) 
one can consider the contribution of the first term in 
the series i)57fl. which comes from the integral involving 
the square of the wave function (|65[) . For the range of 
energies corresponding to the first ten excited states, and 
for rj — 100 we obtain that the first term accounts for 
more than 0.9986 of the whole sum in Eq. (|57|l . 

Fig. [7| shows the exact wave function of the first ex- 
cited state for a = ±oo and rj — 100. The wave function 
was evaluated from expansions i(56"|l and i(59"|) . We ob- 
serve that the wave function is strongly elongated in the 
axial direction and that the region of small r gives rather 
small contribution to the total wave function. In the 
transverse direction the wave function exhibits the expo- 
nential behavior predicted by Eq. (|65fl . These matters 
are illustrated in more details in Figs. 0D and showing, 
respectively, the axial and the transverse profiles of the 
wave function. The figures compare the exact profiles 
with the quasi-one-dimcnsional prediction of Eq. (|65H . 
The approximate curves fit very well the exact wave func- 
tion, except for the transverse profile calculated for z = 0. 
In this case, the small p behavior of the wave function, 
presented in more detail in the inset of Fig. [S] is domi- 
nated by the divergent term l/(27rr), which is absent in 
the approximation l|65|) . 
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FIG. 8: The axial profiles of the wave function presented in 
Fig-Ill evaluated for p = 0, p = 0.08 and p — 0.16. The exact 
profiles (solid lines) are compared with the predictions of the 
quasi-one-dimensional approximation given by Eq. 1651 . All 
lengths are scaled in units of a z = -J h/(pMJ z )- 



C. Quasi-2D regime 

In this section we analyze the properties of the wave 
functions in the regime 77 <C 1. Let us first focus on the 
wave function of the ground state. In this case we apply 
the integral representation (|19f) . where for \£\ ~ 1 and 
77 ^C 1 the main contribution to the integration comes 
from t ~ 1. Expanding the integrand for small rj and 
expressing the result in the physical units we obtain 
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We note that the ground-state wave function depends 
exclusively on the trapping frequency in the z direction, 
thus its size is given approximately by d = \Jh/{poj) 
As for the quasi-one-dimensional traps, we have also 
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FIG. 9: The radial profiles of the wave function presented in 
Fig. |7| evaluated for z = 0, z = 0.5 and z = 1.0. The exact 
profiles (solid lines) are compared with the predictions of the 
quasi-one-dimensional approximation given by Eq. 165H . The 
inset shows the details of z — profile for small values of p. 
All lengths are scaled in units of a z = wh/JJuJz). 



found another representation for the wave function of the 
ground-state. It can be applied, for instance, to derive 
approximate axial and radial profiles, which agree very 
well with the exact wave functions, as we have shown 
in |jj. Expanding l|58|l for small 77 and performing an 



integration over t we arrive at 
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where Kq(x) is a modified Bessel function. 

Fig. ^3 presents the exact ground-state wave function 
for a — ±oo and 77 = 0.01, evaluated from Eqs. I|56l) and 
(|59|l . We note that the ground-state wave function for 
small r is nearly isotropic, while for larger r it is slightly 
elongated in the direction of weaker confinement. The 
anisotropy of a bound-state for 77 <C 1 seems to be larger 
than in quasi-one-dimcnsional traps. 

Let us investigate now the properties of the wave func- 
tions for excited states. The approximate form of the 
wave function can be found from the expansion l|59l) . In 
the regime of energies corresponding to the lowest excited 
states, the sum in Eq. (|59|l is dominated by the first term. 
This results from the asymptotic properties of the con- 
fluent hypergeometric function U(c, 1, x), which decays 
faster in x for larger values of the parameter c. Hence, 
the approximate wave function of excited states in quasi- 
two-dimensional traps is given by 
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Approximate result l|o8fl cannot be directly used for 
p = 0, where the confluent hypergeometric function U 
exhibits logarithmic divergence. Due to the same reason, 
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FIG. 10: Exact wave function r^(r) for two atoms inter- 
acting via s-wave pseudopotential and trapped in a harmonic 
potential with r\ — u±/cu z — 0.01. The figure shows the 
ground state for the scatteri ng length a = ±00. All lengths 
are scaled in units of a z = v/ft/(/zWz). 



the latter approximation is also not valid when r — > 0, 
where the exact wave function behaves as l/(27rr). Simi- 
larly as for the quasi-one-dimensional wave function, the 
leading part of the wave function is located outside the 
region of small p and r, which makes the approximation 
IjHSjl quite accurate. 

To estimate the quality of the approximation (|68J) we 
considered the series (|61|) . where the first term comes 
from the integral of square modulus of l|68() . Within the 
range of energies corresponding to the first ten excited 
states, and for 77 = 0.01, we obtain that the first term of 
(|61|l contributes to more than 0.9985 of the total sum. 

The exact wave function of the first excited state for 
a = ±00 and 77 = 0.01 is presented in Fig. ^J It is 
evaluated from the expansions (|B(t|) and (|59J) . We ob- 
serve that the wave function of the excited state is highly 
anisotropic, and elongated in the radial direction, reflect- 
ing the geometry of the trap. More detailed behavior can 
be deduced from Figs. IT2*) and 1131 showing respectively 
the transverse and the axial profiles of the wave function. 
These figures compare the exact profiles with the quasi- 
two-dimensional approximation given by Eq. (|68|l . In the 
case of transverse profiles all the approximate curves are 
almost indistinguishable from the exact ones, except in 
the region of small p when the approximation l|68|l ceases 
to be valid due to the logarithmic divergence of U. On 
the other hand, the approximate axial profiles for p = 5 
and p — 10 fit very well the exact wave function, while 
for p = 0.1 the approximation (|68|l clearly deviates from 
the exact result. The latter discrepancy can be again at- 
tributed to the logarithmic divergence of (|68|l at small p. 
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FIG. 11: Exact wave function r'I'(r) for two atoms inter- 
acting via s-wave pseudopotential and trapped in a harmonic 
potential with r\ = lu±/uj z — 0.01 The figure presents the first 
excited state for the scattering length a — ±oo. All lengths 
are scaled in units of a z = ^JhJJJuJJT). 
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FIG. 12: The radial profiles of the wave function presented 
in Fig. 1111 evaluated for z — 0, z — 1 and z — 2. The 
exact profiles (solid lines) are compared with predictions of 
the quasi- two-dimensional approximation 16811 . All lengths 
are scaled in units of a z = ^h/(pui z ). 



V. FESHBACH RESONANCES 

In this section we extend our theory to the system of 
two interacting atoms in the presence of Feshbach res- 
onances. The latter technique is widely used in recent 
experiments on ultracold atoms and allows to tune the 
value of the scattering length by changing the strength of 
applied magnetic field. Close to the Feshbach resonance, 
the scattering length can take a very large value or can 
change its sign. Obviously, when a becomes too large, 
the interaction between two atoms cannot be described 
in terms of the standard pseudopotential J3J). To see this 
we recall that (J3J) is obtained from the more general form 
of the pseudopotential [9J 



2ttH 2 tan5 (k) d 
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FIG. 13: The axial profiles of the wave function presented 
in Fig. If fl evaluated for p = 0.1, p = 5 and p — 10. The 
exact profiles (solid lines) are compared with predictions of 
the quasi-two-dimensional approximation 1681 . All lengths 
are scaled in units of a z = \fh/{pu z ). 



in the limit when tan<5o ~ —ka. As it can be easily ver- 
ified, the latter assumption is fulfilled for ka <§; l/(fci?o) 
where Rq is the effective range, whereas application of 
the zero-range pseudopotential requires kRo <C 1 . 

Following the work of Bolda et al. |13|. we introduce 
an effective scattering length, which is defined by 



CZeff(-E) = 



tan<5g(fc) 



(70) 



with k related to the kinetic energy by E = h 2 k 2 /(2/i). 
In general, the effective scattering length is energy depen- 
dent, however, for sufficiently small k, a c s(E) reduces to 
the standard scattering length. Now we can reformulate 
our theory, replacing the standard scattering length a by 
o-cs(E). In this way Eq. (|21|l determining the eigenener- 
gies of the system is substituted by 
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E — E 



(71) 



By solving the latter equation in a self-consistent way 
we obtain the eigenenergies, which are valid for arbitrar- 
ily large values of the scattering length. Moreover, by 
performing the analytic continuation of 1)70(1 to negative 
energies (imaginary k), one can properly account for the 
whole spectrum of bound states [l(| • 

We stress that the validity of the discussed model is 
based on the assumption that the effective range of the 
physical potential is much smaller than the harmonic os- 
cillator length d. In this regime it is justified to use the 
pseudopotential ijfj9"|) . which was derived for free space. 
This assumption also guarantees that at distances com- 
parable to the range of the physical potential, the kinetic 
energy, which enters the pseudopotential (|69fl through k, 
is equal to the total energy E. 
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To apply the concept of the effective scattering length 
to Feshbach resonances, we have to specify the depen- 
dence of the s-wave phase shift on k. The theory 
[25L I26J [23 predicts the following dependence for So 



Obg 



arctan 



7 



E — En 



(72) 



where Sb g is the background phase shift, E m denotes the 
energy of the resonance, and A m is the energy shift due to 
the coupling between the open and closed channels. The 
parameter 7 is a "reduced width" of the resonance which 
is related to the usual width T by T = 2^k |2fj. Close to 
the Feshbach resonance we can assume that the energy 
E m varies linearly with the magnetic field strength 



E m {B) — E m (B — B Tes ), 



where 



E„ — 



dE„ 



dB 



(73) 



(74) 



B reB 



and -B ra s denotes the magnetic field strength at which 
the energy of the closed channel crosses the dissociation 
threshold of the open channel. By combining Eqs. H7U|) 
and H72I) . after some straightforward algebra, we obtain 
the following result for the energy-dependent scattering 
length 



a e «(E) = a bg 



AB 



1 - 



( 1 + t) 



B-IB + E/E' m -AB% 



(75) 

where Ef, = ft 2 / (ma^ ) is the binding energy correspond- 
ing to the background scattering length. The resonance 
width AB and the resonance position Bq are related with 
to the previous parameters by 



and 



Bq = B r€ 



AB = 



A m 
EL 



a hK E' 



(76) 



(77) 



For sufficiently small energies (E <C E{,,E' m Bo) the ef- 
fective scattering length becomes independent of E, and 
Eq. (|75(l reduces to the well-known formula 



a c s(E) = ftbg 



1 



AB 



B-Bn 



(78) 



Employing Eqs. 171|) and l|75|l we have calculated the 
energy spectrum for two 87 Rb atoms interacting close 
to the Feshbach resonance at 100 mT 30]. Fig. ITU 
presents energy levels versus magnetic field in the quasi- 
one-dimensional trap with r\ — 100 and uj z = 5 kHz. 
In addition it shows the values of the magnetic field for 
which the effective scattering length diverges. From the 
plot it is clear that in the considered range of energies the 
position of the resonance changes with the energy. A sim- 
ilar calculation performed for the quasi-two-dimensional 
trap with r\ = 0.01 and io z — 500 kHz is shown on Fig. 1151 
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FIG. 14: Energy spectrum for two 87 Rb atoms versus mag- 
netic field B near the Feshbach resonance at 100 mT. The 
atoms are confined in an axially symmetric trap with ui z = 
5 kHz and lu± — 500 kHz. The dotted line shows the value of 
the magnetic field at which the energy-dependent scattering 
length a c g(E) diverges. 
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FIG. 15: Energy spectrum for two Rb atoms versus mag- 
netic field B near a Feshbach resonance at 100 mT. The atoms 
are confined in an axially symmetric trap with ui z — 500 kHz 
and u± = 5 kHz. The dotted line shows the value of the mag- 
netic field at which the energy-dependent scattering length 
a c ff(E) diverges. 



VI. CONCLUSION 

In summary, we have presented a detailed analysis of 
the system of two interacting atoms confined in an axially 
symmetric harmonic trap. We discussed in detail differ- 
ent regimes, in particular the quasi-one- and the quasi- 
two-dimensional geometries. We have shown that these 
two regimes can be applied already for the traps with 
77 > 10 and T) < 0.1, respectively. In this way we demon- 
strate that, at least on the level of two-atom physics, 
realization of low-dimensional systems does not require 
extremely large (or small) ratios of the transverse to axial 
trapping frequencies. 
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We have applied our analytical results to studying the 
system of two atoms with interaction modified by a Fes- 
hbach resonance. To this end we have utilized the con- 
cept of an energy-dependent scattering length [13}, and 
employing well-known results from the theory of Fesh- 
bach resonances we have derived an explicit formula de- 
termining the energy spectrum in terms of the standard 
parameters describing the resonance. Our results can 
be directly implemented to calculate exact dynamics of 
the two ultracold atoms in harmonic traps with arbitrary 
large trapping frequencies and in the presence of Fesh- 
bach resonances. This is particularly important in the 
context of implementation of quantum information pro- 
cessing in systems of trapped ultracold atoms. 



directly from the definition of the Beta function 

B(x,-1/2) = -2(x-^\b(x,1/2). (A3) 

Inserting the series expansion for B(x, —1/2) into the in- 
tegral in Eq. (|A1|) . and performing the integration term 
by term, we arrive at the final result 



$0) =2 - ln(l + x) 

(2k)\ 



« 



ti ^^ 



(*• 



i)ln- 



(A4) 



Now we turn to the derivation of Eq. (|49|l . We begin 
with the multiplication formula for the function ip{ z ) lit! 
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APPENDIX A: DETAILS OF DERIVATION FOR 
QUASI-TWO-DIMENSIONAL SPECTRUM 

In this appendix we present the derivation of Eqs. I)44|) 
and H49|) , which we used in the analysis of energy spec- 
trum in quasi-two-dimensional traps. Proof of the for- 
mer formula starts from the integral defining the function 
$(x) 



where B(x,y) = r(x)r(y)/r(a; + y) is the Euler Beta 
function. Next, we apply the series representation for 
B(x, 1/2)0 



1 OC 

B(x,l/2) = - + J2 



(2fc)! 1 

x ' ^ (2 k k!) 2 x + k' 



(A2) 



The series expansion of B(x, —1/2) can be found by com- 
bining Eq. (|A2|) with the following identity, which follows 



1 n_1 / k\ 

ip(nz) = - ^ V 1 ( -z + - ) + Inn, (A5) 



fc=0 



and use the definition of B(x,y) to obtain the following 
identity: 

2V?^ rQs + g) , 

> — -, — - — - + w(nx) = 

1 71 — 1 

-£( S (* + 7T>-i)+^+f))+ ln «' 



m=0 



(A6) 



For «>1 we approximate the summation on the r.h.s. 
of Eq. I|A6() by integration. Since the singularities of 
B{x,— 1/2) and ip(x) cancel each other at x = 0, the 
replacement of the summation by integration is valid for 
X > — 1, where the integrand is free from singularities. 
The latter approximation results in 



2, WZ i-^TT + Hnx) 



n>l 



-$(x)-lnx= dt B{x + t,-\) , x>0, (Al) J 

Jo 



n ^ n T(x 

1 

dt(B(x + t,-±)+ip(x + t))+lnn. (A7) 



m=0 * v *" 2 n ■ 

1 



The integral of the first term can be expressed in terms of 
the function <I>(x) (cf. Eq. (|A1|) ~). whereas the integration 
of ip(x) is trivial and follows directly from the definition: 
ip(x) = (d/dx) lnT(x). Finally we obtain 






+ ip(nx) « — O(x) + Inn. 



(A8) 
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